How does urbanization influence ibis behaviors relevant to transmission?

Specifically, how does variation in flock size and frequency of provisioning influence these behaviors?

Data:

10-minute focal follow observations of birds at 5 urban parks in South Florida. Follows conducted during summer 2019.

Focal follows were only conducted in urban parks

Note: Focal follows were only conducted when people were not actively feeding the birds so foraging means probing on the ground for food. The few instances of consuming human provided food are cases where ibis were foraging ‘naturally’ and happened upon scraps of anthropogenic food that was already on the ground.

Sample Sizes / Data Structure

Observations consist of 168 focal follows lasting > 2 minutes. Five urban parks were samples for ~30 focals per site. Sites were visited 2-4 times each.

Provisioning estimates: Each visit, # of people and groups feeding birds were recorded for duration of visit. Provisioning was estimated by counting the number of groups that fed the birds per visit and dividing by total time I was at the site to get # provisioning events/ hour. Each focal follow was paired with the average # of provisioning events per hour on day/site of the focal follow.

Flock size estimates: Every ~30 minutes, flock size was recorded. I paired each focal follow with the nearest estimate of flock size.

Birds seem to spend the most time foraging, grooming, and being vigilant.

I am going to re-group these variables according to their relevance to different modes of parasite transmission (see different sections below).

Summary Statistics

How many focal follows lasted > 2 minutes? *important- when we add in predictors later, we lose data for 3 follows because they lack corresponding predictor data.

## [1] 150

What was the total amount of time spent observing ibis? (this only includes obs >2 minutes aka the ones used in the analysis)

Total obs time in seconds:

## [1] 79994

Total obs time in minutes:

## [1] 1333.233

Mean obs time:

## [1] 533.2933

How many different behaviors were observed on average for each ibis during the observation?

##    focalID           n_behaviors   
##  Length:150         Min.   :2.000  
##  Class :character   1st Qu.:3.000  
##  Mode  :character   Median :4.000  
##                     Mean   :3.853  
##                     3rd Qu.:5.000  
##                     Max.   :6.000

How many birds showed 2 or more behaviors?

nrow(sum3.1[sum3.1$n_behaviors >= 2,]) #all birds observed for at least 2 minutes did 2 or more behaviors 
## [1] 150
nrow(sum3.1[sum3.1$n_behaviors < 2,])
## [1] 0

Most commonly observed behaviors observed?

How many birds foraged?

sum4 <- simp2

sum4f <- subset(sum4, Activity=="Foraging" & activity_count > 0)
summary(sum4f)
##    focalID          focalduration     Activity         activity_count 
##  Length:110         Min.   :139.0   Length:110         Min.   :  4.0  
##  Class :character   1st Qu.:592.0   Class :character   1st Qu.: 53.5  
##  Mode  :character   Median :596.0   Mode  :character   Median :148.0  
##                     Mean   :537.1                      Mean   :212.6  
##                     3rd Qu.:597.0                      3rd Qu.:362.0  
##                     Max.   :598.0                      Max.   :585.0

How many birds groomed?

sum4g <- subset(sum4, Activity=="Grooming" & activity_count > 0)
summary(sum4g)
##    focalID          focalduration     Activity         activity_count 
##  Length:113         Min.   :122.0   Length:113         Min.   :  2.0  
##  Class :character   1st Qu.:557.0   Class :character   1st Qu.: 22.0  
##  Mode  :character   Median :596.0   Mode  :character   Median :140.0  
##                     Mean   :536.5                      Mean   :170.8  
##                     3rd Qu.:597.0                      3rd Qu.:275.0  
##                     Max.   :598.0                      Max.   :572.0

How many birds were vigilant?

sum4v <- subset(sum4, Activity=="Vigilant" & activity_count > 0)
summary(sum4v)
##    focalID          focalduration     Activity         activity_count  
##  Length:140         Min.   :122.0   Length:140         Min.   :  4.00  
##  Class :character   1st Qu.:579.8   Class :character   1st Qu.: 41.75  
##  Mode  :character   Median :596.0   Mode  :character   Median :132.50  
##                     Mean   :539.2                      Mean   :178.65  
##                     3rd Qu.:597.0                      3rd Qu.:269.00  
##                     Max.   :598.0                      Max.   :591.00

What behaviors did ibis spend the most time doing?

sum5<-simp2

sum5 <- sum5 %>% group_by(Activity) %>% summarise(sum(activity_count))
sum5
## # A tibble: 7 × 2
##   Activity `sum(activity_count)`
##   <chr>                    <dbl>
## 1 Drinking                   243
## 2 Foraging                 23387
## 3 Grooming                 19305
## 4 Other                      812
## 5 Resting                   2934
## 6 Vigilant                 25011
## 7 Walking                   8302

Which behaviors were observed most frequently in ibis?

sum6<-simp2

sum6<- as.data.frame(subset(simp2, activity_count > 0))

sum6.1 <- sum6 %>% group_by(Activity) %>% 
  summarise(n_focals=n(),
            .groups = 'drop')

sum6.1
## # A tibble: 7 × 2
##   Activity n_focals
##   <chr>       <int>
## 1 Drinking       19
## 2 Foraging      110
## 3 Grooming      113
## 4 Other          61
## 5 Resting        20
## 6 Vigilant      140
## 7 Walking       115

Visualizations of predictors

Predictors:

What is the distribution of provisioning frequency? What about average flock size? Do we see any relationship between # feeding events and average flock size?

Note: Here I’m showing avg flock size and nearest flock size

Model visualizations and Model Selection

Water Exposure

This variable includes foraging in water, bathing, and drinking.

Main Q’s:

Does flock size influence behaviors related to the exposure of waterborne pathogens?

Does the amount of human provided food influence behaviors related to the exposure of waterborne pathogens?

Exploratory Visualizations

First, what does the distribution of water exposure look like? Plot water exposure by flock size, provisioning, and site to look for potential trends. (set to ‘show nothing, run code’ so rmd isn’t super long)

Models

Standardize predictor variables

Model selection:

Build a super model and then remove variables in a biologically relevant way to see how fit changes and compare AIC values.

Hypothesis testing approach:

Is flock size important? Is provisioning important?

Model diagnostics

Not yet finished. Code is set to ‘don’t run’ at the moment.

Ectoparasite removal

This variable includes preening and bathing.

(set to ‘show nothing, run code’ so rmd isn’t super long)

Foraging / Preening trade-off

Main Q’s:

Does flock size, provisioning, both, or neither influence foraging time? What about preening time?

Literature suggests potential trade-offs between these behaviors. Would be interesting to see if we see differences in either of them and if so, are the effects opposite each other?

Models

Standardize predictor variables

#library(mosaic)  # standardizing variables

ForagePreen_wide <- 
  ForagePreen_wide %>% 
  mutate(avg_flock_size_s = scale(avg_flock_size))%>% 
  mutate(flocksize_nearest_s = scale(flocksize_nearest))%>%
  mutate(feedingperhour_s = scale(feedingperhour))

Foraging models:

Foraging includes foraging on soil, in water, on human objects, and consuming human provided food.

Note: Focal follows were only conducted when people were not actively feeding the birds so foraging means probing on the ground for food. The few instances of consuming human provided food are cases where ibis were foraging ‘naturally’ and happened upon scraps of anthropogenic food that was already on the ground.

Exploratory Visualizations

Model selection:

Build a super model and then remove variables in a biologically relevant way to see how fit changes and compare AIC values. (set to ‘show nothing, run code’ so rmd isn’t super long)

Simple version:
#Simple model selection:

##No random effects: 

#intercept only
summary(fmod0.3)
##  Family: nbinom2  ( log )
## Formula:          Foraging ~ 1 + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1571.0   1580.1   -782.5   1565.0      147 
## 
## 
## Dispersion parameter for nbinom2 family (): 0.987 
## 
## Conditional model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -0.91841    0.09627   -9.54   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.0340     0.1881  -5.497 3.85e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning
fmod11 <-glmmTMB(Foraging~feedingperhour_s+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(fmod11)
##  Family: nbinom2  ( log )
## Formula:          Foraging ~ feedingperhour_s + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1567.9   1579.9   -779.9   1559.9      146 
## 
## 
## Dispersion parameter for nbinom2 family (): 1.03 
## 
## Conditional model:
##                  Estimate Std. Error z value Pr(>|z|)    
## (Intercept)      -0.95398    0.09454 -10.091   <2e-16 ***
## feedingperhour_s -0.19827    0.08581  -2.311   0.0209 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.0324     0.1878  -5.497 3.86e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#flock size
fmod12 <-glmmTMB(Foraging~flocksize_nearest_s+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(fmod12)
##  Family: nbinom2  ( log )
## Formula:          Foraging ~ flocksize_nearest_s + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1567.6   1579.6   -779.8   1559.6      146 
## 
## 
## Dispersion parameter for nbinom2 family (): 1.03 
## 
## Conditional model:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)         -0.94968    0.09438  -10.06   <2e-16 ***
## flocksize_nearest_s -0.26320    0.10240   -2.57   0.0102 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.0343     0.1881    -5.5  3.8e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning + flock size
fmod13 <-glmmTMB(Foraging~flocksize_nearest_s+feedingperhour_s+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(fmod13)
##  Family: nbinom2  ( log )
## Formula:          
## Foraging ~ flocksize_nearest_s + feedingperhour_s + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1565.2   1580.3   -777.6   1555.2      145 
## 
## 
## Dispersion parameter for nbinom2 family (): 1.06 
## 
## Conditional model:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)         -0.97975    0.09296 -10.539   <2e-16 ***
## flocksize_nearest_s -0.24333    0.10301  -2.362   0.0182 *  
## feedingperhour_s    -0.18472    0.08670  -2.131   0.0331 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.0330     0.1878  -5.499 3.82e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning x flock size 
fmod14 <-glmmTMB(Foraging~flocksize_nearest_s*feedingperhour_s+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(fmod14)
##  Family: nbinom2  ( log )
## Formula:          
## Foraging ~ flocksize_nearest_s * feedingperhour_s + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1567.1   1585.2   -777.6   1555.1      144 
## 
## 
## Dispersion parameter for nbinom2 family (): 1.07 
## 
## Conditional model:
##                                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                          -0.97711    0.09345 -10.456   <2e-16 ***
## flocksize_nearest_s                  -0.24949    0.10501  -2.376   0.0175 *  
## feedingperhour_s                     -0.20032    0.09952  -2.013   0.0441 *  
## flocksize_nearest_s:feedingperhour_s -0.05511    0.17401  -0.317   0.7515    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.0333     0.1879  -5.499 3.81e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
AICtab(fmod0.3,fmod11,fmod12,fmod13,fmod14)
##         dAIC df
## fmod13  0.0  5 
## fmod14  1.9  6 
## fmod12  2.4  4 
## fmod11  2.7  4 
## fmod0.3 5.8  3
##Same models but with site as random effect 
#intercept only
summary(fmod0.2)
##  Family: nbinom2  ( log )
## Formula:          Foraging ~ 1 + (1 | site) + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1570.3   1582.4   -781.2   1562.3      146 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance Std.Dev.
##  site   (Intercept) 0.08102  0.2846  
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 1.05 
## 
## Conditional model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -0.9674     0.1599  -6.052 1.43e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)   -1.033      0.188  -5.498 3.84e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning
fmod15 <-glmmTMB(Foraging~feedingperhour_s+(1|site)+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(fmod15)
##  Family: nbinom2  ( log )
## Formula:          
## Foraging ~ feedingperhour_s + (1 | site) + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1569.8   1584.9   -779.9   1559.8      145 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance Std.Dev.
##  site   (Intercept) 0.01202  0.1097  
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 1.04 
## 
## Conditional model:
##                  Estimate Std. Error z value Pr(>|z|)    
## (Intercept)       -0.9640     0.1126  -8.561   <2e-16 ***
## feedingperhour_s  -0.1861     0.1062  -1.752   0.0798 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.0327     0.1879  -5.497 3.86e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#flock size
fmod16 <-glmmTMB(Foraging~flocksize_nearest_s+(1|site)+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(fmod16)
##  Family: nbinom2  ( log )
## Formula:          
## Foraging ~ flocksize_nearest_s + (1 | site) + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1569.6   1584.6   -779.8   1559.6      145 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance  Std.Dev. 
##  site   (Intercept) 4.094e-07 0.0006398
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 1.03 
## 
## Conditional model:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)         -0.94967    0.09439 -10.062   <2e-16 ***
## flocksize_nearest_s -0.26319    0.10244  -2.569   0.0102 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.0343     0.1881    -5.5  3.8e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning + flock size
fmod17 <-glmmTMB(Foraging~flocksize_nearest_s+feedingperhour_s+(1|site)+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(fmod17)
##  Family: nbinom2  ( log )
## Formula:          
## Foraging ~ flocksize_nearest_s + feedingperhour_s + (1 | site) +  
##     offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1567.2   1585.3   -777.6   1555.2      144 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance  Std.Dev. 
##  site   (Intercept) 2.283e-09 4.779e-05
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 1.06 
## 
## Conditional model:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)         -0.97975    0.09296 -10.539   <2e-16 ***
## flocksize_nearest_s -0.24333    0.10301  -2.362   0.0182 *  
## feedingperhour_s    -0.18472    0.08670  -2.131   0.0331 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.0330     0.1878  -5.499 3.82e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning x flock size 
fmod18 <-glmmTMB(Foraging~flocksize_nearest_s*feedingperhour_s+(1|site)+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(fmod18)
##  Family: nbinom2  ( log )
## Formula:          
## Foraging ~ flocksize_nearest_s * feedingperhour_s + (1 | site) +  
##     offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1569.1   1590.2   -777.6   1555.1      143 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance  Std.Dev. 
##  site   (Intercept) 2.117e-09 4.601e-05
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 1.07 
## 
## Conditional model:
##                                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                          -0.97710    0.09345 -10.455   <2e-16 ***
## flocksize_nearest_s                  -0.24949    0.10501  -2.376   0.0175 *  
## feedingperhour_s                     -0.20032    0.09952  -2.013   0.0441 *  
## flocksize_nearest_s:feedingperhour_s -0.05511    0.17401  -0.317   0.7515    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.0333     0.1879  -5.499 3.81e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
AICtab(fmod0.2,fmod15,fmod16,fmod17,fmod18) #fmod17 and fmod18 = best fit
##         dAIC df
## fmod17  0.0  6 
## fmod18  1.9  7 
## fmod16  2.4  5 
## fmod15  2.6  5 
## fmod0.2 3.1  4
Pseudo R^{2} calculation

I think you could do this using the MuMIn package. It calculates conditional and marginal coefficient of determination for Generalized mixed-effect models. Marginal represents the variance explained by the fixed effects and conditional is interpreted as a variance explained by the entire model including both fixed and random effects.

library(MuMIn)
## 
## Attaching package: 'MuMIn'
## The following object is masked from 'package:bbmle':
## 
##     AICc
r.squaredGLMM(fmod17,fmod0.2) #warning: effects of zero-inflation and dispersion model are ignored 
## Warning: 'r.squaredGLMM' now calculates a revised statistic. See the help page.
## Warning in r.squaredGLMM.glmmTMB(fmod17, fmod0.2): the effects of zero-inflation
## and dispersion model are ignored
##                   R2m         R2c
## delta     0.028590314 0.028590314
## lognormal 0.063815282 0.063815284
## trigamma  0.007691387 0.007691387

Ben Bolker also wrote a work-around for doing this directly. This can be accessed using the ‘performance’ add-on package (or using his crude code on github)

Link to discussion for future reference: https://github.com/glmmTMB/glmmTMB/issues/169

##THESE METHODS SOMEHOW DON'T ACCOUNT FOR RANDOM EFFECTS
performance::r2(fmod18)
## Warning: Can't compute random effect variances. Some variance components equal
##   zero. Your model may suffer from singularity (see `?lme4::isSingular`
##   and `?performance::check_singularity`).
##   Solution: Respecify random structure! You may also decrease the
##   `tolerance` level to enforce the calculation of random effect variances.
## Warning: mu of 0.4 is too close to zero, estimate of random effect variances may
##   be unreliable.
## Random effect variances not available. Returned R2 does not account for random effects.
## # R2 for Mixed Models
## 
##   Conditional R2: NA
##      Marginal R2: 0.008
performance::r2_nakagawa(fmod17)
## Warning: Can't compute random effect variances. Some variance components equal
##   zero. Your model may suffer from singularity (see `?lme4::isSingular`
##   and `?performance::check_singularity`).
##   Solution: Respecify random structure! You may also decrease the
##   `tolerance` level to enforce the calculation of random effect variances.

## Warning: mu of 0.4 is too close to zero, estimate of random effect variances may
##   be unreliable.
## Random effect variances not available. Returned R2 does not account for random effects.
## # R2 for Mixed Models
## 
##   Conditional R2: NA
##      Marginal R2: 0.008
performance::model_performance(fmod18)
## Random effect variances not available. Returned R2 does not account for random effects.
## # Indices of model performance
## 
## AIC      |     AICc |      BIC | R2 (cond.) | R2 (marg.) |    RMSE | Sigma | Score_log | Score_spherical
## --------------------------------------------------------------------------------------------------------
## 1569.121 | 1569.909 | 1590.195 |            |      0.008 | 171.286 | 1.065 |   -11.037 |           0.042
performance::model_performance(fmod12)
## Random effect variances not available. Returned R2 does not account for random effects.
## # Indices of model performance
## 
## AIC      |     AICc |      BIC | R2 (cond.) | R2 (marg.) |    RMSE | Sigma
## --------------------------------------------------------------------------
## 1567.594 | 1567.870 | 1579.636 |            |      0.005 | 177.225 | 1.029
performance::r2(fmod13)
## Warning: mu of 0.4 is too close to zero, estimate of random effect variances may
##   be unreliable.
## Random effect variances not available. Returned R2 does not account for random effects.
## # R2 for Mixed Models
## 
##   Conditional R2: NA
##      Marginal R2: 0.008
Hypothesis testing approach:

Is flock size important? Is provisioning important? (How important are they?) (set to ‘show nothing, run code’ so rmd isn’t super long)

Model diagnostics

Testing the fit of fmod17 and fmod18- the two best fitting models.

First test the fit of fmod17

##  Family: nbinom2  ( log )
## Formula:          
## Foraging ~ flocksize_nearest_s + feedingperhour_s + (1 | site) +  
##     offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1567.2   1585.3   -777.6   1555.2      144 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance  Std.Dev. 
##  site   (Intercept) 2.283e-09 4.779e-05
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 1.06 
## 
## Conditional model:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)         -0.97975    0.09296 -10.539   <2e-16 ***
## flocksize_nearest_s -0.24333    0.10301  -2.362   0.0182 *  
## feedingperhour_s    -0.18472    0.08670  -2.131   0.0331 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.0330     0.1878  -5.499 3.82e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

## 
##  DHARMa nonparametric dispersion test via sd of residuals fitted vs.
##  simulated
## 
## data:  simulationOutput
## dispersion = 0.65756, p-value = 0.148
## alternative hypothesis: two.sided

## 
##  DHARMa zero-inflation test via comparison to expected zeros with
##  simulation under H0 = fitted model
## 
## data:  simulationOutput
## ratioObsSim = 1.0019, p-value = 1
## alternative hypothesis: two.sided

Residuals look fine. Plotting residuals against the two predictors suggests that there might be some issues for one of them.

Now test fit for fmod18

##  Family: nbinom2  ( log )
## Formula:          
## Foraging ~ flocksize_nearest_s * feedingperhour_s + (1 | site) +  
##     offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1569.1   1590.2   -777.6   1555.1      143 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance  Std.Dev. 
##  site   (Intercept) 2.117e-09 4.601e-05
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 1.07 
## 
## Conditional model:
##                                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                          -0.97710    0.09345 -10.455   <2e-16 ***
## flocksize_nearest_s                  -0.24949    0.10501  -2.376   0.0175 *  
## feedingperhour_s                     -0.20032    0.09952  -2.013   0.0441 *  
## flocksize_nearest_s:feedingperhour_s -0.05511    0.17401  -0.317   0.7515    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.0333     0.1879  -5.499 3.81e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

## 
##  DHARMa nonparametric dispersion test via sd of residuals fitted vs.
##  simulated
## 
## data:  simulationOutput
## dispersion = 0.66608, p-value = 0.182
## alternative hypothesis: two.sided

## 
##  DHARMa zero-inflation test via comparison to expected zeros with
##  simulation under H0 = fitted model
## 
## data:  simulationOutput
## ratioObsSim = 1.0021, p-value = 1
## alternative hypothesis: two.sided

The residuals look pretty okay for this model too but there are more problems when I plot against specific predictors.

Preening models:

Preening includes preening or grooming.

Exploratory Visualizations

Model selection:

Build a super model and then remove variables in a biologically relevant way to see how fit changes and compare AIC values. (set to ‘show nothing, run code’ so rmd isn’t super long)

Simple version:
#Simple model selection:

##No random effects: 

#intercept only
summary(pmod0.3)
##  Family: nbinom2  ( log )
## Formula:          Preening ~ 1 + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1545.9   1554.9   -769.9   1539.9      147 
## 
## 
## Dispersion parameter for nbinom2 family (): 0.763 
## 
## Conditional model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)   -1.158      0.108  -10.72   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.1912     0.2037  -5.849 4.95e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning
pmod11 <-glmmTMB(Preening~feedingperhour_s+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(pmod11)
##  Family: nbinom2  ( log )
## Formula:          Preening ~ feedingperhour_s + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1547.2   1559.3   -769.6   1539.2      146 
## 
## 
## Dispersion parameter for nbinom2 family (): 0.766 
## 
## Conditional model:
##                  Estimate Std. Error z value Pr(>|z|)    
## (Intercept)      -1.16616    0.10803 -10.795   <2e-16 ***
## feedingperhour_s  0.09753    0.12266   0.795    0.427    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.1919     0.2038  -5.848 4.98e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#flock size
pmod12 <-glmmTMB(Preening~flocksize_nearest_s+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(pmod12)
##  Family: nbinom2  ( log )
## Formula:          Preening ~ flocksize_nearest_s + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1545.0   1557.0   -768.5   1537.0      146 
## 
## 
## Dispersion parameter for nbinom2 family (): 0.783 
## 
## Conditional model:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)          -1.1765     0.1067 -11.022   <2e-16 ***
## flocksize_nearest_s   0.1638     0.1016   1.613    0.107    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.1870     0.2027  -5.856 4.74e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning + flock size
pmod13 <-glmmTMB(Preening~flocksize_nearest_s+feedingperhour_s+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(pmod13)
##  Family: nbinom2  ( log )
## Formula:          
## Preening ~ flocksize_nearest_s + feedingperhour_s + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1546.6   1561.6   -768.3   1536.6      145 
## 
## 
## Dispersion parameter for nbinom2 family (): 0.784 
## 
## Conditional model:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)         -1.18239    0.10677 -11.074   <2e-16 ***
## flocksize_nearest_s  0.15642    0.10063   1.554    0.120    
## feedingperhour_s     0.07797    0.11909   0.655    0.513    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.1875     0.2028  -5.855 4.76e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning x flock size 
pmod14 <-glmmTMB(Preening~flocksize_nearest_s*feedingperhour_s+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(pmod14)
##  Family: nbinom2  ( log )
## Formula:          
## Preening ~ flocksize_nearest_s * feedingperhour_s + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1548.6   1566.6   -768.3   1536.6      144 
## 
## 
## Dispersion parameter for nbinom2 family (): 0.784 
## 
## Conditional model:
##                                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                          -1.18402    0.10912 -10.851   <2e-16 ***
## flocksize_nearest_s                   0.15785    0.10248   1.540    0.123    
## feedingperhour_s                      0.08129    0.12786   0.636    0.525    
## flocksize_nearest_s:feedingperhour_s  0.01284    0.17902   0.072    0.943    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.1877     0.2029  -5.855 4.77e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
AICtab(pmod0.3,pmod11,pmod12,pmod13,pmod14)
##         dAIC df
## pmod12  0.0  4 
## pmod0.3 0.9  3 
## pmod13  1.6  5 
## pmod11  2.2  4 
## pmod14  3.6  6
##Same models but with site as random effect 
#intercept only
summary(pmod0.2)
##  Family: nbinom2  ( log )
## Formula:          Preening ~ 1 + (1 | site) + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1545.8   1557.8   -768.9   1537.8      146 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance Std.Dev.
##  site   (Intercept) 0.07115  0.2667  
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 0.802 
## 
## Conditional model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.2169     0.1624  -7.491 6.84e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.1857     0.2025  -5.856 4.73e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning
pmod15 <-glmmTMB(Preening~feedingperhour_s+(1|site)+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(pmod15)
##  Family: nbinom2  ( log )
## Formula:          
## Preening ~ feedingperhour_s + (1 | site) + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1547.7   1562.8   -768.9   1537.7      145 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance Std.Dev.
##  site   (Intercept) 0.08273  0.2876  
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 0.805 
## 
## Conditional model:
##                  Estimate Std. Error z value Pr(>|z|)    
## (Intercept)      -1.22063    0.17031  -7.167 7.65e-13 ***
## feedingperhour_s -0.03225    0.17501  -0.184    0.854    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.1851     0.2023  -5.857 4.72e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#flock size
pmod16 <-glmmTMB(Preening~flocksize_nearest_s+(1|site)+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(pmod16)
##  Family: nbinom2  ( log )
## Formula:          
## Preening ~ flocksize_nearest_s + (1 | site) + offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1546.6   1561.7   -768.3   1536.6      145 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance Std.Dev.
##  site   (Intercept) 0.03272  0.1809  
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 0.799 
## 
## Conditional model:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)          -1.2026     0.1384  -8.686   <2e-16 ***
## flocksize_nearest_s   0.1329     0.1219   1.091    0.275    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.1848     0.2022  -5.858 4.68e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning + flock size
pmod17 <-glmmTMB(Preening~flocksize_nearest_s+feedingperhour_s+(1|site)+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(pmod17)
##  Family: nbinom2  ( log )
## Formula:          
## Preening ~ flocksize_nearest_s + feedingperhour_s + (1 | site) +  
##     offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1548.5   1566.6   -768.2   1536.5      144 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance Std.Dev.
##  site   (Intercept) 0.01502  0.1226  
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 0.793 
## 
## Conditional model:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)         -1.19393    0.12748  -9.366   <2e-16 ***
## flocksize_nearest_s  0.14585    0.11664   1.250    0.211    
## feedingperhour_s     0.05801    0.15247   0.380    0.704    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.1861     0.2026  -5.855 4.77e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning x flock size 
pmod18 <-glmmTMB(Preening~flocksize_nearest_s*feedingperhour_s+(1|site)+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=ForagePreen_wide)
summary(pmod18)
##  Family: nbinom2  ( log )
## Formula:          
## Preening ~ flocksize_nearest_s * feedingperhour_s + (1 | site) +  
##     offset(log(duration))
## Zero inflation:            ~1
## Data: ForagePreen_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1550.5   1571.6   -768.2   1536.5      143 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance Std.Dev.
##  site   (Intercept) 0.01476  0.1215  
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 0.792 
## 
## Conditional model:
##                                       Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                          -1.194963   0.128629  -9.290   <2e-16 ***
## flocksize_nearest_s                   0.147157   0.118857   1.238    0.216    
## feedingperhour_s                      0.060848   0.161153   0.378    0.706    
## flocksize_nearest_s:feedingperhour_s  0.009679   0.179544   0.054    0.957    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.1862     0.2026  -5.855 4.78e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
AICtab(pmod0.2,pmod15,pmod16,pmod17,pmod18)
##         dAIC df
## pmod0.2 0.0  4 
## pmod16  0.8  5 
## pmod15  2.0  5 
## pmod17  2.7  6 
## pmod18  4.7  7

The top performing models are within 2 AIC values of the intercept only null model.

Hypothesis testing approach:

Is flock size important? Is provisioning important? (How important are they?) (set to ‘show nothing, run code’ so rmd isn’t super long)

Model diagnostics

Not needed because the intercept only model is the best fitting model?

Vigilance models:

#library(mosaic)  # standardizing variables

Vig_wide <- 
  Vig_wide %>% 
  mutate(avg_flock_size_s = scale(avg_flock_size))%>% 
  mutate(flocksize_nearest_s = scale(flocksize_nearest))%>%
  mutate(feedingperhour_s = scale(feedingperhour))

Exploratory Visualizations

Model selection:

Build a super model and then remove variables in a biologically relevant way to see how fit changes and compare AIC values. (set to ‘show nothing, run code’ so rmd isn’t super long)

Simple version:
#Simple model selection:

##No random effects: 

#intercept only
summary(vmod0.3)
##  Family: nbinom2  ( log )
## Formula:          Vigilant ~ 1 + offset(log(duration))
## Zero inflation:            ~1
## Data: Vig_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1795.1   1804.1   -894.5   1789.1      147 
## 
## 
## Dispersion parameter for nbinom2 family (): 1.04 
## 
## Conditional model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.1186     0.0831  -13.46   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -2.7527     0.3674  -7.493 6.72e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning
vmod11 <-glmmTMB(Vigilant~feedingperhour_s+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=Vig_wide)
summary(vmod11)
##  Family: nbinom2  ( log )
## Formula:          Vigilant ~ feedingperhour_s + offset(log(duration))
## Zero inflation:            ~1
## Data: Vig_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1795.2   1807.3   -893.6   1787.2      146 
## 
## 
## Dispersion parameter for nbinom2 family (): 1.05 
## 
## Conditional model:
##                  Estimate Std. Error z value Pr(>|z|)    
## (Intercept)      -1.12883    0.08280 -13.634   <2e-16 ***
## feedingperhour_s  0.11559    0.08451   1.368    0.171    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -2.7571     0.3691  -7.471 7.98e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#flock size
vmod12 <-glmmTMB(Vigilant~flocksize_nearest_s+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=Vig_wide)
summary(vmod12)
##  Family: nbinom2  ( log )
## Formula:          Vigilant ~ flocksize_nearest_s + offset(log(duration))
## Zero inflation:            ~1
## Data: Vig_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1796.7   1808.7   -894.3   1788.7      146 
## 
## 
## Dispersion parameter for nbinom2 family (): 1.05 
## 
## Conditional model:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)         -1.11978    0.08298 -13.494   <2e-16 ***
## flocksize_nearest_s -0.05246    0.08319  -0.631    0.528    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -2.7504     0.3665  -7.504 6.21e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning + flock size
vmod13 <-glmmTMB(Vigilant~flocksize_nearest_s+feedingperhour_s+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=Vig_wide)
summary(vmod13)
##  Family: nbinom2  ( log )
## Formula:          
## Vigilant ~ flocksize_nearest_s + feedingperhour_s + offset(log(duration))
## Zero inflation:            ~1
## Data: Vig_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1796.6   1811.7   -893.3   1786.6      145 
## 
## 
## Dispersion parameter for nbinom2 family (): 1.06 
## 
## Conditional model:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)         -1.13097    0.08263 -13.688   <2e-16 ***
## flocksize_nearest_s -0.06498    0.08119  -0.800    0.423    
## feedingperhour_s     0.12392    0.08560   1.448    0.148    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -2.7546     0.3682  -7.482 7.33e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning x flock size 
vmod14 <-glmmTMB(Vigilant~flocksize_nearest_s*feedingperhour_s+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=Vig_wide)
summary(vmod14)
##  Family: nbinom2  ( log )
## Formula:          
## Vigilant ~ flocksize_nearest_s * feedingperhour_s + offset(log(duration))
## Zero inflation:            ~1
## Data: Vig_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1798.5   1816.6   -893.2   1786.5      144 
## 
## 
## Dispersion parameter for nbinom2 family (): 1.06 
## 
## Conditional model:
##                                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                          -1.13733    0.08423 -13.503   <2e-16 ***
## flocksize_nearest_s                  -0.05786    0.08320  -0.695    0.487    
## feedingperhour_s                      0.13835    0.09467   1.461    0.144    
## flocksize_nearest_s:feedingperhour_s  0.05327    0.14684   0.363    0.717    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -2.7558     0.3687  -7.474 7.75e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
AICtab(vmod0.3,vmod11,vmod12,vmod13,vmod14)
##         dAIC df
## vmod0.3 0.0  3 
## vmod11  0.2  4 
## vmod13  1.5  5 
## vmod12  1.6  4 
## vmod14  3.4  6
##Same models but with site as random effect 
#intercept only
summary(vmod0.2)
##  Family: nbinom2  ( log )
## Formula:          Vigilant ~ 1 + (1 | site) + offset(log(duration))
## Zero inflation:            ~1
## Data: Vig_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1796.7   1808.8   -894.4   1788.7      146 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance Std.Dev.
##  site   (Intercept) 0.01613  0.127   
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 1.06 
## 
## Conditional model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.1253     0.1008  -11.16   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -2.7493     0.3662  -7.508 5.99e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning
vmod15 <-glmmTMB(Vigilant~feedingperhour_s+(1|site)+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=Vig_wide)
summary(vmod15)
##  Family: nbinom2  ( log )
## Formula:          
## Vigilant ~ feedingperhour_s + (1 | site) + offset(log(duration))
## Zero inflation:            ~1
## Data: Vig_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1796.9   1812.0   -893.5   1786.9      145 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance Std.Dev.
##  site   (Intercept) 0.0167   0.1292  
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 1.07 
## 
## Conditional model:
##                  Estimate Std. Error z value Pr(>|z|)    
## (Intercept)      -1.13354    0.10079 -11.247   <2e-16 ***
## feedingperhour_s  0.12838    0.09711   1.322    0.186    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -2.7515     0.3671  -7.495 6.61e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#flock size
vmod16 <-glmmTMB(Vigilant~flocksize_nearest_s+(1|site)+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=Vig_wide)
summary(vmod16)
##  Family: nbinom2  ( log )
## Formula:          
## Vigilant ~ flocksize_nearest_s + (1 | site) + offset(log(duration))
## Zero inflation:            ~1
## Data: Vig_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1798.4   1813.4   -894.2   1788.4      145 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance Std.Dev.
##  site   (Intercept) 0.01607  0.1268  
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 1.06 
## 
## Conditional model:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)         -1.12738    0.10085  -11.18   <2e-16 ***
## flocksize_nearest_s -0.05547    0.09087   -0.61    0.542    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -2.7480     0.3657  -7.514 5.73e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning + flock size
vmod17 <-glmmTMB(Vigilant~flocksize_nearest_s+feedingperhour_s+(1|site)+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=Vig_wide)
summary(vmod17)
##  Family: nbinom2  ( log )
## Formula:          
## Vigilant ~ flocksize_nearest_s + feedingperhour_s + (1 | site) +  
##     offset(log(duration))
## Zero inflation:            ~1
## Data: Vig_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1798.5   1816.6   -893.2   1786.5      144 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance Std.Dev.
##  site   (Intercept) 0.01066  0.1032  
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 1.07 
## 
## Conditional model:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)         -1.13444    0.09477 -11.971   <2e-16 ***
## flocksize_nearest_s -0.05885    0.08814  -0.668    0.504    
## feedingperhour_s     0.12762    0.09276   1.376    0.169    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -2.7515     0.3671  -7.495 6.64e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#provisioning x flock size 
vmod18 <-glmmTMB(Vigilant~flocksize_nearest_s*feedingperhour_s+(1|site)+offset(log(duration)),
                    family=nbinom2,
                    ziformula=~1,
                    data=Vig_wide)
summary(vmod18)
##  Family: nbinom2  ( log )
## Formula:          
## Vigilant ~ flocksize_nearest_s * feedingperhour_s + (1 | site) +  
##     offset(log(duration))
## Zero inflation:            ~1
## Data: Vig_wide
## 
##      AIC      BIC   logLik deviance df.resid 
##   1800.4   1821.5   -893.2   1786.4      143 
## 
## Random effects:
## 
## Conditional model:
##  Groups Name        Variance Std.Dev.
##  site   (Intercept) 0.007085 0.08417 
## Number of obs: 150, groups:  site, 5
## 
## Dispersion parameter for nbinom2 family (): 1.06 
## 
## Conditional model:
##                                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                          -1.13790    0.09207 -12.359   <2e-16 ***
## flocksize_nearest_s                  -0.05544    0.08769  -0.632    0.527    
## feedingperhour_s                      0.13647    0.09856   1.385    0.166    
## flocksize_nearest_s:feedingperhour_s  0.03804    0.16156   0.235    0.814    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Zero-inflation model:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -2.7533     0.3679  -7.483 7.24e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
AICtab(vmod0.2,vmod15,vmod16,vmod17,vmod18)
##         dAIC df
## vmod0.2 0.0  4 
## vmod15  0.2  5 
## vmod16  1.6  5 
## vmod17  1.8  6 
## vmod18  3.7  7

The top performing models are within 2 AIC values of the intercept only null model.

Hypothesis testing approach:

Is flock size important? Is provisioning important? (How important are they?) (set to ‘show nothing, run code’ so rmd isn’t super long)

Model diagnostics

Not needed because the intercept only model is the best fitting model?